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Abstract 

A simple analysis method based upon a transverse shear deformation theory 
and a sublaminate approach is utilized to analyze a Mixed-Mode edge delamination 

specimen. The analysis provides closed form expressions for the interlaminar 

r 

shear stresses ahead of the crack, the total energy release rate, and the erfergy 
release rate components. The parameters controlling the behavior are identffied. 
The effect of specimen stacking sequence and delamination interface on the strain 
energy release rate components is investigated. Results are compared with a 
finite element simulation for reference. The simple nature of the method makes 
it suitable for preliminary design analyses which require a large number of 
configurations to be evaluated quickly and economically. 


* This work was sponsored by the NASA Langley Research Center under Grant 
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Introduction 


Delaminations along the free edges of laminates subjected to tensile 

loading have been observed during testing and service. The presence of 

delamination, initiated by interlaminar stresses, causes redistribution of 

the stresses among plies in a laminate and, therefore, usually results in a 

reduction of stiffness and strength. The edge delamination (ED) test has 

been proposed by Pagano and Pipes 1 to characterize the interlaminar peel 

2 

strength of laminated composite materials. O'Brien extended the scope of 

the test, to investigate delamination onset and growth in graphite/epoxy 

laminates udder uniform extension. . A simple expression was also developed 

for the total energy release rate. The energy release rate components Gj, 

Gjj and G^. associated with the opening, shearing and tearing modes, 

respectively, were estimated based on a finite element simulation and the 

2 3 

crack-closure method ’ . 

A similar approach was used to study delaminations around an open hole 
4 

in composite laminates . Discrete locations around the hole boundary were 
modeled as straight edges, with the ply orientations rotated by an 
appropriate angle. Delamination was found to be governed by the percentage 
of Mode I for a given geometry under static loading for the graphite/epoxy 
material systems under consideration. 

5 

Whitney and Knight developed an ED specimen which produces Mode I 
behavior. The analysis was based on classical laminated plate theory and 

v ^ 

continuity of displacements, force resultants, and moment resultants 
between the cracked and uncracked regions of the plate were not satisfied. 
In addition, such an approach precludes any reasonable determination of the 
effect of specimen geometry on Gj. 


> 


g 

Recently Whitney developed a higher order laminated plate theory which 

includes transverse shear deformation and a thickness-stretch mode to 

analyze a Mode I ED specimen. The effect of specimen geometry on strain 

energy release rate was also investigated for 12 ply laminates of the class 

[e/-e 2 /e/90 2 3 s and for the class [0 3 /90 3 ] s . 

It is the purpose of the present work to develop a simple model for 

the analysis of Mixed-Mode ED specimens. Such a model provides closed-form 

estimates of Gj, and Gjjj and hence allows one to establish appropriate 
7 8 

failure criteria ’ for delamination. 

Preliminary Remarks 

Consider the ED specimen shown in Figure 1 subjected to a uniform 
strain e = e. Due to symmetry one quarter of the laminate is analyzed as 

A 

shown in Figure 2. The response is only a function of Y and Z. The 

laminate is divided into sublaminates with thicknesses h. and h and local 

1 o 

coordinate sys.tems z^, y^ and z q , y Q . The crack length is denoted by a. 

Sublaminates 1 and 2, and 0 and 3 represent the groups of plies above and 

below the interface along which delamination occurs, respectively. 

In order to provide an accurate estimate of interlaminar stresses, a 

higher-order theory should be considered since classical laminated plate 

theory predicts zero interlaminar stresses. A shear deformation theory can 

be used for this purpose. This theory provides a good estimate for 

interlaminar shear stresses x and x . However, the interlaminar peel 

xz yz 

stress o zz is not accurate. The reason for that is the absence of 
thickness strain. 
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From symmetry, the transverse displacement w is zero at 1=0 , hence the 
prescription of > w at the middle plane Z=0 fixes w everywhere. In this 
case, the vertical shearing force resultant at both ends cannot be 
prescribed and the distribution of the peel stress will not be correct. 
Inspite of this simplification, reliable energy release rate components can 
be estimated based on interlaminar shear stresses. Gj is evaluated as - 
(Gii+Gin), where the total energy release rate is denoted by G^. 

In the present formulation thickness strain is neglected and 
consequently considerable simplification in the' analysis is achieved. 

Another source of simplification in the present approach is due to the 

modeling of the structure as sublaminates — group of plies that are 

conveniently treated as laminated units. This approach can be applied with 

confidence if the characteristic length of the response is large compared 
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to the individual sublaminate thickness. 

Overview of the Analytical Solution 

In the following sections a step by step procedure is provided for the 
solution of the ED specimen. Intermediate results are also provided. The 
governing equations are derived in Appendix I. Expression for the 
interlaminar stresses, total energy release rate and energy release rate 
components are given in Equations (34), (40) and (42)-(43), respectively. 
The parameters associated with these equations are provided explicitly in 
terms of the stiffness coefficients in Appendix II. 

A solution based on Classical Lamination Theory is given in Appendix 
III. This solution represents the behavior in the interior of the 
laminate. Application of the present analysis and comparison with a finite 
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element simulation for 63 test cases is presented under the section 
entitled Results and Discussion. 

The reader interested in results and comparison can refer directly to 
the "Results and Discussion" section on page 28. 

Analysis 

Assume the following displacement field within each sublaminate: 

>■. ' 

u = xe + U(y) + zP x (y) 

v = V(y) + zp y (y) (1) 

w = W(y) 

where u,v, and w denote displacements relative to the x, y, and z axes, 

respectively, and e is a uniform axial strain. Coordinates y and z are 

local coordinates as shown in Figure 2. The present formulation recognizes 

shear deformation through the rotations B and B . The corresponding 

x y 

strains are 


xx 



+ Z£ 


y.y 


e 

Y 

Y 

Y 


ZZ 

xy 

XZ 

yz 


= o 

= U,„ + zp 


x,y 


= P, 


p y + w ’y 


( 2 ) 


The variables associated with sublaminates 0 through 3 will be written with 
subscripts 0 through 3, respectively. 
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From symmetry 


or 


"o (y >' ,, o / 2 > = 0 


W = 0 
o 


( 3 ) 


From continuity of displacements at the interface between sublaminate (1) 
and (0) 

u Q (y, h Q /2) = u 1 (y,-h 1 /2) 
v Q (y, h Q /2) = u 1 (y,-h 1 /2) 
w Q (y, h Q /2) = w 1 (y,-h 1 /2) 

(4) 

Substitute from Equation (1) into (4) to get 

* * , , 

h 1 h n 

U o = U 1 " T P lx " ~2 P ox 


V 


0 



h _l *_0 

2 3 lx 2 P ox 


(5) 

Governing Equations 

The governing equations for each sublaminate are derived in Appendix I 
using a virtual work approach. These equations are written below for 
convenience. 
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Constitutive Relationship 
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e 
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xy 

= 
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6 11 

B 12 

B 16 
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°16 
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B 22 
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°26 


P x,y 

M 

xy J 


B 16 

B 26 

B 66 
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°66_ 



i — — i 
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— 



1 

•O 

<< 

1 


A 44 

A 45 


B y + W ’y 

L°xJ 


_ A 45 
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Equilibrium Equations 


N + n = 0 
xy,y x 

N + n = 0 

y.y y 

Q + q = 0 

y,y 


xy.y 


Q + m = 0 
x 


M - Q + m = 0 

y.y y y 


( 6 ) 


(7) 


where n , n , q, m and m are defined in Equation (1-9) of Appendix I. 
x y x y 

The equilibrium equations can be written in terms of kinematic 
variables by substitution from Equation (6) into Equation (7). The result 
is 
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A, C L 
66 yy 





U ' 


n 

X 

A oc L 
26 yy 

A 22 L yy 


SYM 


V 


n 

y 

0 

0 

~ A 44 L yy 



w 

= - 

-q 

B 66 L yy 

B 26 L yy 

_A 45 L y 

{D 66 L yy " A 55 ) 




m x 

B 26 L yy 

B 22 L yy 

“ A 44 L y 

(D 26 L yy " A 45 ) (D 22 L yy ' A 44 } 


p y 


m y 


L = d 2 /dy 2 

yy J 

L y = d/dy 


( 8 ) 


Solution Methodology 

The above collection of equations are to be applied' to individual 
plies of a laminate or to groups of plies... sublaminates. The solution 
steps are summarized in the following: 

1. Divide the laminate into sublaminates according to geometry and 
loading condition. The sublaminate length is selected such that 
within the sublaminate the geometry and loading are continuous as 
is commonly done in engineering analysis of simple structures. 

2. The displacements, resultant forces and moments, and interlaminar 
stresses in each sublaminate are governed by the equilibrium 
Equation(7), the constitutive relations (6) and the 
displacement distributions (1). Write these equations for each 
sublaminate in the analysis model. 

3. Apply interlaminar continuity conditions and enforce traction or 
displacement conditions at the extreme upper and lower surfaces 
of the laminate. 


8 



4 


4. Solve the system of coupled ordinary differential equations for 
the element variables. 

5. Enforce the boundary conditions at constant values of y, the 

* 8 

laminate^sections, as well as continuity requirements between 
sublaminate ends in order to find the values of the arbitrary 
constants resulting from the solution in step 4. 

6. Determine interlaminar stresses, resultant forces and moments 
displacement distribution and energy release rate. 

Application to the ED Specimen 

The ED configuration is divided into- four sublaminates as shown in 
Figure 2. The response associated with sublaminates 1 and 0 is coupled 
through the continuity conditions at their common interface. Hence, the 
variables associated with both sublaminates are to be solved 

simulataneously. The situation is different with sublaminates 2 and 3 

where the continuity conditions are relaxed due to the presence of the 

crack. Therefore, the variables associated with these sublaminates are not 
coupled. 

The solution procedure for sublaminates 1 and 0--the uncracked portion 
of the laminate — is presented first, followed by the cracked portion 

represented by sublaminates 2 and 3. 

Uncracked Region of the Laminate: 

(i) Sublaminate 1 

The upper surface of this sublaminate is stress free. Denote the 

shear and peel stress at the bottom surface by t^, t and p, respectively, 
The transverse displacement W^, is zero from Equation (5). Hence, the 

equilibrium equation in terms of the displacement variable takes the form 
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J., 

l 66 yy 

A oc L 
26 yy 

B cc L 
66 yy 

B 26 L yy 

i 0C 1 L 
26 yy 

A 

22 L yy 

B 0C 1 L 
26 . yy 

B 22 ll yy 

0 

0 

" A 45 ll y 

“ A 44 lL y 

; 1 L 
66 L 

B 26 ll yy 


)(D 26 lL yy- A - 

^6 ll yy 

B 2 2 \y 

(D 26 1| -yy" A 45 1 

)(D 2 2 V A < 


u 

V 

p 

p 


1 

1 

lx 

iy 


(9) 



(ii ) Sublaminate 0 

From symmetry condition at the sublaminate bottom surface the shear 

stresses t^ x and are zero. From reciprocity of stresses at the interface 

between sublaminates 0 and 1, the interlaminar stresses at the upper surface 

of sublaminate 0 are it , t and p. Denote the peel stress at the sublaminate 

x y 

bottom surface by p^. Hence the equlibrium equation takes the form 


A,,°L 
66 yy 

A ' °L 
26 yy 

B 66°*"yy 

B 26° L yy 


U o 


t 

X 

A °L 
26 L yy 

A °L 
22 L yy 

B 26° L yy 

B 22° L yy 


V o 


l y 

0 

0 

' A 45° L y 

- A 44° L y 
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p OX 

= - 

Pj-P 

B 66° L yy 

B 26° L yy 

( D 66 0L yy- fl 55 0 > 

< D 26° L yy- A 45°» 


8 

p oy 


K 

B 26° L yy 

B 22° L yy 

< D 26 0 V A 45 0 ) 

< D 22° L yy- A 44 0 ) 


.K. 


( 10 ) 
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Equations (9),. (10) and (5) can be combined to yield the following 

form: 


(F nVV) 

F 21 L yy 

^ F 31 L yy" A 45 1) 

F /11 L 
41 yy 


(F 22 L yy" A 44 ) 
F 32 L yy 

(F 42 L yy' A 45 0) 


SYM 


p iy 

(F 33 L yy _A 55 \ 


^oy 


Plx- 

F 43 L yy (F 44 L yy' A 55 } 


^ox 


= 0 


( 11 ) 

The parameters in Equation (11) are defined in Equation ( 1 1 -1 ) of Appendix 
II. 


Assume an exponential solution of the form 


( p ly> p 0 y’ p lx’ P ox^ = ^ly ’ P oy ’ P lx ’ 




( 12 ) 


Substitute from Equation (12) into Equation (11) to get the following 
characteristic equation 

EgS 8 + E 6 s 6 + E 4 s 4 + E 2 s 2 + E q = 0 (13) 

Coefficients E g through Eq are defined in Equation ( 1 1 -5 ) of Appendix II. 

Parameter E Q represents the determinant of the coefficients of L in matrix 

o yy 

equation (11) while Eq is the determinant when is set to zero. Parameter Eq 

depends solely on the stiffness coefficients A^ 4 , A gg and A^ g for both 

sublaminates while E g is predominantly influenced by the bending and coupling 

coefficients D-. and B... Hence its numerical value can be orders of magnitude 
1 J » J 

smaller than the remaining coefficients. This fact results in the presence of a 
boundary zone in the response. 
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Typical values of a nondimensional form of the coefficients Eg through 
Eg is given in Table II for three laminates made of T300/5208 
graphite/epoxy material. The nondimensional form of the coefficients is 
obtained by making the following substitution 


s= sb 


(14) 


Where b is the laminate semi-width. The arrows in the layups of Table II 
indicate ..the interfaces containing delaminations. The material properties 
and geometry N appear in Table I. 

The characteristic roots controlling the behavior are determined from 
Equation (13), which has a closed-form solution. However when A^g, B^g, 
(j = 1,2) and A^g are neglected. Equation (11) takes an uncoupled form 
and consequently the characteristic equation (13) can be factorized into 
two biquadratic equations. The stiffness coefficients A^g, B^g and D^g 
represent coupling interaction while A^g depends predominantly on the shear 
modulii and G^g. For a composite material where these shear modulii 
are approximately the same, A^g can be neglected. 

The uncoupled form of the characteristic equation is 


^ F 11 F 22 _F 21. _ ^ F 11 A 44 + F 22 A 44 + A 44 A 44 ^ 


^ F 33 F 44 _F 43 ~( F 33 A 55 + F 44 A 55 + A 55 A 55 ^ = 0 


( 15 ) 
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f/ The first bracket in Equation (15) control V and behavior while the 
second U and B . The absolute values of the roots in Equation (15) will be 

A 

denoted by s^ , s^. ar >d S X 2 ‘ These can be regarded as a good 

approximation for the roots of the coupled equation (13). A comparison 
between the coupled and uncoupled roots for a typical laminate is provided 
in Table III. These roots are found to be real for the material system and 
layup used. The predictions of the uncoupled equation (15) are in good 
approximation with the coupled equation (13). Also, the larger 
characteristic roots s^ and s ^ correspond to bending behavior while the 
smaller s^ and s^ control the membrane behavior. This is shown in the 
following section. 

Membrane behavior can be modeled by setting 


M 



(16) 


in the equilibrium equations for sublaminates 1 and 0. The characteristic 
equation (13) reduces to 


(F llm s " A 44 )(F 


22m s2 ' A 55 } 0 


(17) 


The membrane parameters F^ and are defined in Equation (II -9 ) 
of Appendix II. They depend on A.j and B.^ coefficients. The 
characteristic roots predicted by Equation (17) are included in Table III. 
By comparison with the roots of Equation (15), bending behavior is more 
localized than the membrane behavior as the characteristic roots 
controlling be.nding are larger. This fact is expected since Classical 
Lamination Theory (CLT), which predicts membrane-type behavior, prevails 
in the interior of the laminate. 
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Since the laminate width is large compared to its thickness and to the 
crack length, the response in sublaminates 1 and 0 is predominantly 
decaying from the crack tip and, therefore, only the roots with a negative 
sign will be considered in this solution. 

Cracked Region of the Laminate: 

Sublaminate 2 

This sublaminate represents the upper group of plies in the cracked 
portion of the laminate. Since there is no restriction on the transverse 
displacement W, boundary conditions on Q can be specified. 

The upper and lower surfaces in this sublaminate are stress-free and 
at-' y^-a there is a' free edge. The equilibrium equation (7) reduces to 




2,y 




(18) 


By substituting these conditions into the constitutive relations (6), to a 
single differential equation in terms of P 2x obtained. 


(D 66 + B 26 lcd 12 + B 66 1 Cd 22 + D Z6' C(i 3Z ) P 2x,yy 


A 


(A V 
1 lrt 45 1 

55 “ a 1 

"44 -» 


*2x = 0 


(19) 

Parameters Cd^, Cd^ and Cd^ are defined in Equation ( 1 1 -10 ) of Appendix 
II. Solution of Equation (18) leads to 


where 


H Vl - s c y l 

P 2x ~ H i e H 2 e 


0 < y. < - a 


14 


s 

c 


A 


55 


( A 45 ^ ^ A 44 


1/2 


^ee 1 + B 26 lcd 12 + B 66 lcd 22 + 


( 20 ) 


Since the crack length 'a' can be small, positive as’ well as negative 
signs of the root s c have been considered. The arbitrary constants and 
H 2 are determined from the boundary condition at 


"'h " ' a M xy 2 <- a > ' 0 < 21 > 

and continuity conditions at y^ = 0 with sublaminate 1. 

The displacements at the bottom surface of sublaminate 2 are 


h l h l 

v 2 (y l’ " h l /2 * = V i^°^ + ^ Cd ll " 2 Cd 31^ ey l ‘ ^ Cd 12" 2 Cd 32 ) 


s c y l 


"Vl, 


[H 1 (1-e A ) + H 2 (l-e " A )] 

h l s c y l ' s c y l 

u 2 (y 1 -h 1 /2) = u.(0) + Cd 21 eyi - (Cd 22 - -f) [Hjd-e ") + H 2 (l-e ")] 


(22) 

where 

* v.(0) = v 2 (0,Ch 1 /2) 

u- (0) = u 2 (0 - h x /2) 


The linear terms in y 1 in Equation (22) represent the displacement 
when p 2x = 0. 
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(ii) Sublaminate 3 

The upper surface of this sublaminate is stress-free, while the 
interlaminar shear stresses on the lower surface are zero from symmetry 
conditions. Moreover, there is a free edge at y Q = -a. The equilibrium 
equation (7) reduces to 



M 

xy 3 ,y 




(23) 


By substituting these conditions into the constitutive relations (6), 
two coupled differential equations in and p 3x is obtained 


L - a 0 
22 yy rt 44 

J 26 L yy" A 45° 

L - A 0 
26 yy H 45 

J 66 L yy" A 55_ 




= 0 

(24) 


Parameters J 22 , J 2 g and J g g 
(11-12) of Appendix II . 

Solve the differential 


in Equation (24) are defined in Equation 


equations for and p 3x 


, s Vo 


3y 


+ I, 


■ ss i y o 


+ I. 


, SS 2 y o 


+ 1/ 


SS 2 y o 


3x 


= n- 


(I. 


ss iV _ 

e + I, 


e ~ SS l y O) 


d 3 e 


ss 2 y o 


~ ss 2 y o 

I 4 e ) 


where 


-a < y 


o 


< 0 


16 


(25) 


"J * ' < J 22 Ss j ‘O' ‘ J 26 ss j - A 45°> 


The constants 1^ through 1^ are found from the boundary conditions at 


Y = -a (-a) = M (-a) = 0 (26) 

o y - 3 x y 3 

and continuity conditions at y =0 with sublaminate 0. 

o 

The roots ss^ and ss^ in Equation (25) are found by solving the 
characteristic equation resulting from equation (24). 

The displacements at the upper surface of sublaminate 3 are given by 


h ss.y 

v 3 (y o > h Q /2) = v.(o) + wd 12 ey Q + (wd 12 + j- + HjWd^) [I^e ■ °- 1) 

Ti 


+ I 2 (e 


' SS l y o 


SS 2 y o 


1)] + (wd 12 + ~2 + r ) 2 wd 13 ) CI 3 ( e " !) 


~ ss 2 y o 

+ I 4 (e • - D] 


h ss iy n 

u 3 (y 0 , ho/2) = u^o) ■+ wd 21 ey 0 +' [wd 22 + ( wd 23 + g 2 ) ] [^(e - 1) 

" ss l y 0 h o SS 2 y 0 

+ I 2 (e - 1)] + [wd 22 + n 2 ( wd 23 + _ 2 )] CI 3 (e - 1) 

~ ss 2 y o 

+ I 4 (e * 0 - 1)] (27) 

Parameters wd^, wd^, wd^, wd^, wd 22> and wd^ are defined in Equation 
(11-12) of Appendix II . 

In order to determine the energy release rate components by the 
virtual crack-closure method 1 ^, the relative displacements at the crack 
surface as well as the interlaminar stresses at the crack tip are needed. 
From Equations (22) and (27) the relative displacements are 
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Av = v 3 C^o- h o /2) ' v 2 (y l ’ " h l /2) 

h. h 

= (wd^ - Cd 11 + 2 “ cd 3i) + ( wd i2 + “2 + n l wd i3 ) 

ss.y -ss.y 

[I x (e - 1) + I 2 (e - 1)] 

b ss 0 y -ss 0 y 

+ (wd 12 + ^ + n 2 wd i3^ n 3 (e - 1) + I 4 (e - 1)] 

b 1 s y -s y 

+ (Sd 12 - ^ Cd 32 ) [H 1 (1 - e c ) + H 2 (1 - e c )] 

au ^ ,u 3 (y 0 , h Q /2) - u 2 (y x , - hj/2) 

\ 

h o ss l y 

= (wd 21 - Cd 21 )ye + [wd 22 + ( wd 23 + 2“)] C I 2 ( e ~ 1 ) 

-ss.y 

+ I 2 (e " D] 

h ss„y -ss„y 

+ ^ wd 22 + n 2 ^ wd 23 + 2 ~^ ^3 " 1) + I4 ( e _ 1)] 

h, s y -s y 

+ (Cd 22 - jt) [Hjd - e c ) + H 2 (1 - e C )] 

-a < y < 0 
(28) 

The linear terms in Equation (28) represent the relative displacements 
when the shear deformations g 2x , g 3x and p 3y are neglected. The remaining 
terms are exponential, and their effect oh the predictions are depicted in 
Figures 4 and 5, for two typical laminates. The dotted lines denoted by 
AU^ and AV^ in the figures correspond to the linear contribution, while the 
solid lines AU and A V include the exponential terms. The crack length 'a 1 
in figures is 10 percent of the laminate semi-width and the applied strain 
is 1000 micro in/in; The simple linear distribution represents a good 
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approximation for the test cases considered. Due to the simplicity and 
accuracy of the linear displacement distribution, the energy release rate 
computations will be based on the linear displacement contribution only. A 
discussion of the effect of the exponential terms on the calculation of the 
total energy release rate will be provided later. 

Interlaminar Stresses 

Denote the absolute values of the roots in Equation (13) by s. (j= 

J 

1,4). The response can be written in the form 




y 


B. = G -a.e 
lx J J 


-s,y 


B = G . e 
p oy rj 


-s,y 


B = G ■ y • e 
ox J T J 


-s-y 


-s .y 

U. =. C- - s.v.G.e J 
l.y ■ v j j j 


V y = C U - s j u j G j e 


■s .y 


j=l,4;o<y<b 

(29) 


Summation over the range of index j is implied in Equation (29). 

Parameters a., y., v., and u. are defined in Equation (11-14) of 

J J J J 3 

Appendix II. The arbitrary constants of integration associated with the 
displacements are C and C^. 
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The parameters associated with Equation (30) are defined in Equation 
(11-15) of Appendix II . 

There are twelve arbitrary constants of integration namely: G, (j = 

J 

1,4), C v , C u , H^, and 1^ through 1^. These can be found from the 

boundary condition at the free edge of sublaminates 2 and 3 . expressed in 
Equations (21) and (26) and from the • continuity conditions between 
sublaminates 1 and 2, and 0 and 3 at y = 0. These are: 




N (0) = N (0) = 0 


N , (0) = N (0) = 0 


xy 


l 


xy c 


M (0) = M (0) = 0 
y l y 2 

N„ (0) = N (0) = 0 


N . (0) 

= N ( 

X 

*< 

o 

xy 3 

M , (0) 

= M ( 

xy 2 v 

xy x v 

M , (0) 

= Ku ( 

xy o 

xy 3 

M (0) = 

= l\ (0) 

y 

^0 

y 3 


*2x<°> = Pi x <°) 


P 0 x<°> = P 3 x<°> 

P 0 y(°) = V° } 


(31) 
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Boundary .conditions on cannot be specified in sublaminates 1 and 0 
as a result of neglecting the transverse normal strain in the assumed 
displacement function. Equations (21), (26) and (31) represent 14 boundary 
conditions for twelve constants. However, the last two continuity 
conditions at the interface between sublaminates 0 and 3 in Equation (31) 
lead to 


that is, . 


44 


V°> 


+ A 


45 


0 = A 44° V°) + A 


45 


e3x<°> 


Qy <°> = Qy (0) 
*0 y 3 


(32) 


which cannot be specified. Therefore, the conditions P QX (0) = P 3x (0) ancl 
3 Qy (0) .= (0) cannot be prescribed for consistency. 

Substitute for the resultant forces and moments in sublaminates 0 
through 3 into the first nine equations in (31) and solve for the arbitrary 
constants. As and are easily obtained, their expressions are listed 
below for convenience. 


where 


(k 26 k 16 
(k 26 k 12 


k 66 k 12 ^ £// ° 
k 22 k 16^ e/ ^ D 



= A, ^ 

+ A. ~ 

12 

"12 

12 


= A. c ^ 

+ A, 

16 

"16 

16 


(33) 


Parameters kgg, 

(II-4). 


k^ and 0 are defined in Equations ( 1 1 -3 ) and 
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The resultant axial forces and moments in Equation (30) consist of two 
parts: a constant term and an exponentially decaying term. In the 
interior of the laminate the second term is negligible and the response is 
controlled by the first term. Since CLT is retrieved in the interior of 
the laminate, the first term can be recognized as the CLT prediction. This 
is shown in Appendix III where a CLT solution is derived for the ED 
specimen. This'' solution is simple and can be derived on a ply-by-ply 
basis. Moreover, the sign of the resultant axial force N y determines the 
sign of the interlaminar peel stress. A compressive peel stress tend to 
retard delamination at a given interface. 

Interlaminar stresses at the interface between sublaminates 1 and 0 
are found from equilibrium. 


t = N. 


2 “ S i y 
N s. G.e J 


x xy l ,y x lj J J 

2 ' s i y 

t = N = N . s. G. e J 
Y Y x ,y YjJ J J 


p = Q y p y ' s j (A 44 1 + a 45 J a j } G j e Sjy 


(34) 


The distribution of the peel stress p is not in equilibrium since the 
boundary condition on the shear force Q,(0) cannot be prescribed. 




pdy = 


Qu „ dy = Q w (b) - Q (0) 


0 y i’ y 


r l 


(35) 


t 0 . 


Q (b) = 0 as Q is an exponentially decaying function, however Q ( 0 ) 
y l y y l 
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A comparison of the interlaminar shear stress distribution t along 
the interface between sublaminates 1 and 0 appears in Figure 6 for three 
laminates. The applied uniform strain is 1000 micro in/in. 

Energy Release Rate 


The total energy release rate can be determined by considering the 

5 

work done by external forces . For a uniform applied strain e x = e 
throughout the laminate this can be written as 


\ 


r _ -e dP 

g t “ T dl 


(36) 


where P is the axial applied force and is given by 

,(b-a) r o 

P = (N - + N ) dy + (N + N )dy 

o X 1 o J -a x 2 x 3 


(37) 


The axial forces N and N in sublaminate 1 and 0 respectively are 

X 1 x o 

given in Equation (30) while the axial forces in sublaminate 2 and 3 are 
given by 


N 


x2 


11 1 1 s c y -V 

(A 11 + A 12 Cd ll + A 16 Cd 21 + B 12 Cd 31 )e + J l S c (H l e ' H 2 e } 


ss 2 y -ss : y 


x3 


(A 11 + A 12 wd ll + A 16 wd 21 )e + (J 2 + n l J 3 )SS 1 (I 1® l Z e } 


ss 2 y -ss 2 y 
+ (J 2 + n 2 J 3 ) ss 2 (I 3 e - I^e ) 

(38) 
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where 


= S 16 1 

+ A ^ 
H 12 

Cd 12 + A 16 Cd 22 + B 12 Cd 32 


= B 12° 

+ A 0 
rt 12 

wd 12 + A 16 wd 22 


= B 16° 

+ a 0 
"13 

wd 13 + A 16 wd 23 

(39) 


Substitute from Equations (30) and (38) into Equation (36) and using 
Equation (37) to get 


G T ~ 2 [K 12 C v + K 16 C u " * A 12* Cd ll + '' * 16 * Cd 21 + B 12 * Cd 31 )E 


' 1 


- (A 12 ° wd + A 16 ° wd 21 ) e] + G r 


(40) 


where G R is an exponential function of the crack length. For a crack 
length larger than a few ply thicknesses, the contribution of G R is 
negligible. This situation is depicted in Figure 7 where a normalized 
total energy release rate for a typical laminate is plotted against the 

r 

crack length divided by ply thickness. The laminate thickness is denoted 
by t in the figure and the ply thickness is h. In this case, Gy reaches 
the constant value predicted by the first term in Equation (40) for a value 
of crack length larger than four ply thickness. Also appearing in Figure 7 
is a comparison with a finite element solution presented in Figure 6 of 
Reference 4. 

The total energy release rate is a global parameter which does not 
depend strongly on the local details at the crack tip. This is the reason 
why relatively simple modeling approaches yield adequate predictions for 
the total energy release rate. 
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This finding checks the results of References 4 and 8 where simple 

> 

closed-form expressions for Gy are in good agreement with predictions based 
on refined theories or on finite element simulation. 


Energy Release Rate Components 


Using the virtual crack-closure method and Gjjj are given by 


G rT = Lim 

1 1 n ' rt Cv) 

. -v 6 h>-0 


,6 

t (S-r)Av(r)dr 
0 y 


1 f 5 

G ut = Lim t (6-r)Au(r)dr 

111 6+0 ■ J 0 


(41) 

where 6 is a virtual crack step size. Unless a singularly exists in the 
stress field, Equation (41) yields the trivial result (Gjj= Gjjj = 0) when 

9 

the limit as <5 tends to zero is determined. Consequently, a sufficiently 
large finite step size is essential to get an answer when using models that 
do not exhibit singular behavior. 

Subsitute from Equation (34) for the interlaminar stresses in term of 
force resultants into Equation (41) and use the linear terms in Equation 
(28) for the relative displacements to get 


G II “ 2 0 v F y 


G I 1 1 = 2 0 u F x 


26 


where 


0 

0 


v 

u 



wd n“ cd n + 2 cd 3i 


wd 21 - Cd 21 


N w dy 


0 y l 



■ 0 



dy 

1 


(42) 


The finite crack step size is denoted by A. The resultant forces N . 
and N ^ exhibit a boundary layer behavior. This is shown in Figures 8 and 
9 where the resultant force distributions are plotted along the laminate 
width. In the interior a constant value corresponding to the CLT 
prediction is reached. 

The average resultant force values F and F depends on the selected 

y * 

value of the crack step A. Recommended ad hoc values as a percentage of 
the initial crack length have been suggested in the finite element 
representation of the crack-closure method. 1 ^ However, within the boundary 
layer region N ^ and N ^ have steep gradients and consequently a small 
variation in the selection of A leads to large variations in F and F . 

y * 

Since and distributions at the crack tip are controlled by the 

boundary layer decay length, the crack step size A should be selected based 
on the boundary layer length rather than a percentage of the initial crack 
1 ength. 

Mode I energy release rate is found from 

Gj = G-p - (Gjj + Gjjj) (43) 
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Gjj and Gj.jj are given in Equation (42) and G^. in Equation (40). 

The energy release rate components predicted by this approach show a 
good correlation with finite element solutions when the crack step size A 
is selected as 




* + 


(44) 


where s. ( j =1 , 4 ) are the characteristic roots in increasing, magnitude 
J 

nondimensional ized by the ply thickness h. The distinct roots s. control 

* J 

the decay length associated with different physical variables. It has been 
found empirically that the combination above gives very good correlation 
with finite element simulations for over sixty cases that have been 
compared. Although the functional form is not simply identified with a 
boundary length it nevertheless contains the proper information. 

For a general layup the crack size A is influenced by all four roots. 
However, for laminates where Mode III is negligible the crack step size is 
influenced by the characteristic roots s 2 and s^ that control V and 
behavior. For this situation, the following crack step size expression 
provides good correlation with the finite element solution: 


I = 2.6 (s 4 i 2 ) I 2 (45) 
Results and Discussion 

An extensive comparison between the energy release rate components 
predicted by the present approach and a quasi 3-D finite element solution 
has been performed. The results appears in Tables IV through X and in 
Figures 10 through 28. The finite element results in Tables V through X 



are related to. the work of Reference 4 and were provided by the authors. 
The properties and geometry used appear in Table I. The applied uniform 
strain is one micro in/in. 

In Table IV a comparison is provided for six laminates. The first 
three layups were reported in Reference 3. The Gjjj component in these 
layups are negligible. The remaining layups represent cases where Gj, G^ 
and Gjjj components are all finite. 

The laminates presented in Tables V-X are rotated stacking sequences 
for two quasi-isotropic layups. The first is a [0/90/+45] s and the second 
is a [45/90/-45/0] s . The results of the first laminate appear in Tables 
V-VII while in Tables VIII-X the results of the second laminate are 
provided. These laminates were used in Reference 4 to investigate 
delaminations around an open hole. The strain energy release rate 
distribution around the hole boundary, for delaminations growing in a 
prescribed interface, was calculated by assuming that each circumferencial 
portion acts as a straight edge subjected to an appropriate uniform 
circumferencial strain. Hence, at each circumferencial angle 0 the 
laminate reflects a new stacking sequence where the load is applied in the 
0-direction tangent to the hole. The angular position corresponding to 
each stacking sequence is provided in Tables V-X. 

The results in Table V are plotted in Figures 10-12. In Figures 10 
and 11 the percentages of Gjj and Gjjj are plotted against the angular 
position, 0, around a hole in the first laminate. The finite element 
results are shown in solid lines while the results of the present analysis 
appear in dotted lines. A comparison of G^ is shown in Figure 12. 
Similarly, the results in Table VI and VII are plotted in Figures 13-16 and 
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17-19, respectively. The results of the second laminate are* provided in 
Tables VIII-X. The results in Table VIII are shown in Figures 20-22. 
Similarly the results of Tables IX and X are plotted in Figures 23-25 and 
26-28. In Figures 22, 25 and 28 as well as in Tables VIII-X the values of 
G-p from the simple expression derived by O'Brien are included for 
comparison. The plots in Figures 26-28 have been discontinued at 0 = -45° 
as compressive peel stress occurs at the 45/90 interface. The sign- of the 

O ‘ 

V 

peel stress at a given interface can be determined from a simple 
membrane-type model. This is shown in Appendix III. Under compressive 
peel stress v the crack surfaces tend to close and a special modeling 
approach should be used. One possible approach was proposed in Reference 9 
in connection with a double cracked-lap-shear specimen tested under 
compression loading. 

Conclusion 

Interlaminar stresses and energy release rates are estimated for the 
ED test using a shear-type deformation theory and a sublaminate approach. 
The predictions are obtained in closed form and the parameters controlling 
the behavior are identified. The governing equations are derived using a 

i 

virtual work approach. Due to the absence of transverse strain the 
interlaminar peel stress distribution is not in equilibrium. The 
interlaminar shear stresses, however, show reliable trends. The energy 
release rate components Gjj and Gjjj are estimated based on the 
interlaminar shear stresses and relative displacements using the virtual 
crack closure method. The total energy release, Gy, is determined from the 
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rate of change of the work done by the external forces with crack length. 
Then Gj is obtained as the difference between Gy and (Gjj + Gjjj)- 

An extensive comparison between the energy release rates predicted by 
the present approach and a quasi 3-D finite element solution for over sixty 
test cases, is performed. The agreement is good. The methodology outlined 
in this work is simple and the results are generated using a Hewlett 
Packard 9845B desktop computer. 
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Appendix I 

Derivation of the Governing Equations 


In this Appendix the governing equations for the sublaminate shown in Figure 
3 are derived using the principle of virtual work. 

Consider a sublaminate... a single ply or group of plies conveniently 
treated as laminated units of thickness h. The origin of a cartesian coordinate 
system is located within the central plane (x-y) with the z-axis being normal to 
this plane. The material of each ply is assumed to possess a plane of elastic 
symmetry parallel 'to xy as shown in Figure 3. 

Stress and moment resultants, 


(N x , N , N xyJ Q x , Q ) = 

h/2 

I <v v T «- V dz 


h/2 

-h/2 

(M , M , , M ) = 

v x ’ y’ xy' J 

( v 

<j , t )zdz 
y’ xy' 


-h/2 



Because of the existence of a plane of elastic symmetry, the constitutive 
relations are given by 


a 

X 


c n 


£ x 

a y 


C 12 C £2 SYM 


E y 

°z 


C 13 C 23 C 33 


E z 

1 

X 

<< 

1 


_ C 16 C 26 C 36 C 66_ 


Y xy 


T 


C-. SYM 


Y 

yz 

= 

44 


T yz 

N | 
X 

1 


_^45 C 55 


_ Y x^ 


( 1 - 2 ) 
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where C,. are components of the anisotropic stiffness matrix and y , y 

1 J '** J ^ 

and y are engineering shear strains 

A Z 

The displacements are assumed to be of the form 
u = U(x,y) + zp x (y) 
v = V(y) + zp y (y) 

w = W(y) (1-3) 


where u,.v and w are the displacement components in the x, y and z 
directions, v respectively. Equation (1-3) in conjunction with the 
strain-displacement relations of classical theory of elasticity leads to 
the following kinematic relations 


xx 


'yy 


u 'x 

% + zp y.y 


zz 


= 0 


r xy V ’y + zP x,y 


'xz 


■yz 


= K 


p y + W ’y 


v 


d-4) 


Substitute Equation (1-4) into Equation (1-2) and put the results into 
Equation ( I - 1 ) . This yields the following constitutive 
rel ations: 
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where 


(A ij’ B i j ’ D ij } ‘ | C ij (1, z ’ z )dz 

-h/2 


The usual' coupling between bending and extensional modes for a laminate of 

i 

arbitrary construction occurs in Equation (1-5) through the stiffness terms 



The variation of the strain energy due to virtual displacements 6u, 6 v 
and 6w is 
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where Se^, 6e^, 6e z> 6Y X y, <$Yy Z and <S Y XZ are the strains associated with 
the virtual displacements. Substitute from Equation (1-3) and integrate 
through the thickness using Equation (1-1) to get 


60 = [N 60, + N 6V, + N 6U, + Q 63 + Q (63 + 6W, ) + M 63 

J fl L x ’x y y xy ’y y x p x y y ' p y ’y' x p x,; 


+ M 63 + M 6 3 ] dA 

y y,y xy p x,y J 


(1-7) 


The variation of the work done by the external forces and by the 
surface fractions is 


6W 


= (n 6U + n 6V + q6W + m 63 + m 63 ) dA 
j a A y a a y y 


* 

(N 60 + N 60 + M 63 + M 63 )ds 
c n n ns s n n ns s' 

•' X 


( 1 - 8 ) 


where a bar denotes values on the boundary, n and s are coordinates normal 
and tangential to the edge, and 


"x = ‘zx ' ‘lx 


n = try ~ t. 

y 2y ly 
q = P 2 - P x 


m = o (t 0 + t. ) 

x 2 v 2x lx' 


"y = 2 <‘2, + V 


(1-9) 


where n and n can be regarded as effective distributed axial forces, m 
x y x 

and m y effective distributed moments and q an effective lateral pressure. 
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From the principle of virtual work the equations of equilibrium and 
boundary conditions are determined from the Euler equations and boundary 
conditions of the variational equation. 


6 V = 6W 


(1-10) 


Substitution of Equations (1-7) and (1-8) into Equation ( I - 10 ) leads the 
following equations of equilibrium: 


N + N + n = 0 
x,x xy,y x 


N + N + n = 0 
xy,x y,y y 


Q + Q + q = 0 
y x,x v y,y M 


M +M - Q + m = 0 

x,x xy,y v x x 

M + M - Q + m = 0 

xy,x y,y v y y 


( 1 - 11 ) 


and one member of the following five products must be prescribed on the 
sublaminate edges 

N n V N ns U s' Vn- M ns and «n W 

For the ED specimen under uniform extension, U(x,y) in Equation (1-3)- 
is given by 


U(x,y) = xe + U*(y) (1-13) 

and the response is a function of y and z coordinates only. For this case 
the equilibrium equations ( I -11 ) take the form 
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0 


N 

xy,y 


+ n x 


= 0 


y,y 


+ n 


= 0 


Q + q = 0 

w y,y M 


M - Q + m = 0 
xy ,y y x x 


M - Q + m = 0 

y,y v y y 


N 


(1-14) 


Substitution of the constitutive relations in Equation (1-5) into Equation 


A,.L 
66 yy 


A 26 L yy 

A 22 L yy 


0 

0 

“ A 44 L yy 

B k< -L 
66 yy 

B 26 L yy 

“ A 45 L y 

B, C L 
26 yy 

B 22 L yy 

~ A 44 L y 




u* 

SYM 


V 



w 

< D 66 L yy- A 55> 


p x 

<D Z6 L yy" A 45* (D 22 L yy' A 44 ) 


p y 






(1-15) 


where the operators 



Ly = d/dy 


and for the boundary conditions at y = constant prescribe 

N or U*, N or V, Q or W, M or 8 and M or p 
y y y y a j a 


(1-16) 


(1-17) 
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Appendix II 

Definition of Parameters 


In this appendix the parameters in Equations (11) - (30) are defined 
in terms of the sublaminates thickness and A.., B. . and D.. coefficients. 

J J 1 J 

Sublaminates 1 and 0 : 

The F.j parameters ( i =1-3 , j =1-4 ) in Equation (11) are defined as 

V* KS 

F 11 = °22 1 x+ T B 22 + h 22 V 11 + h 26 U 11 
F 21 = h 22 v 12 - h 26 u 12 

F 31 = D 26 + 2~ B 26 + h 22 V 13 + h 26 U 13 

F 41 = h 22 v i4 + h 26 u 14 
h 

F 22 °22 “ T B 22 + C 22 V 12 + C 26 U 12 

F 32 = h 26 v 12 + h 66 u 12 
F 43 = h 26 v 14 + h 66 u l4 

F 33 = D 66 1 + T + h 26 V 13 h 66 U 13 

F 42 = D 26° ' T B 26° + C 22 V 14 + C 26 U 14 

F 44 = °66 T 8 66 + C 26 V 14 + C 66 U 14 (II-l) 
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Parameters u^. and (i= 1-4) in Equation ( 1 1 -1 ) relate the 
displacements and \l ^ to the rotations through the following equation 


v i,yy 


_ u i,yy„ 



V 11 L yy v 12 L yy v 13 L y 

u n1 L... u,„L u 10 L 


v l4 L yy 


p iy 

s oy 

u 14 L yy 


S lx 

_ 


^ox 




where 

V 11 = ^ K 26 h 26 " K 66 h 22^ /D + *V 2 
v 12 = ^ K 26 C 26~ K 66 C 22^ D + *V 2 

V 13 = ,^ K 26 h 66 " K 66 h 26^ /D 
v 14 = (K 26 C 66 " K 66 C 26- /D 
U 11 = ^ K 26 h 22 " K 22 h 26^ /D 

u 12 = ^ K 26 C 22" K 22 C 26^ D 

U 13 = ^ K 26 h 26" *22 h 66^ D + h l' /2 

u 14 = ^ K 26 C 26" K 22^66^ D + h o /2 


(II-2) 


D = K 22 K 66 ' K 26 


( I I“3) 
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and 


h - a 1 + _i A 1 
n 22 °22 2 h 22 


h 26 = 8 26 1 + 2 k 2§ 


1 1 1 
h = R ■*■ + —i. A A 

n 66 “66 2 *66 


r - o i 0 a 3 
L 22 " “22 2 h 22 


h 

r = r 0 + _2. a 1 

°26 “26 2 rt 26 


C 66 = B 66° + T A 66 

K 22 = A 22 2 + A 22° 

K 26 = A 26 + A 26 

K 66 = A 66 + A 66 


(II“4) 


The coefficients Eg through in Equation (13) are defined as 


Eg = F n W 1 - F 31 x 1 + F 21 Yj - F 41 Z x 

E 6 = F 11 W 2 " A 44 W 1 “ F 31 X 2 + A 45 X X 1 + F 21 Y 2 " F 41 Z 2 
E a = CA 44 ° A 55 ° - (A., 0 ) 2 ] (F„F„ - F, 2 ) (F„ A? K + F dA A° 


45 


11 33 31' 


22 55 44 44 


2F 42 A 45° } (F 11 A 55 ' F 31 A 45 1} " A 44* W 2 + X 2 + F 21 Y 3 

E o = '^ A 44 A 55 " ( A 45 ) 2 1 t F 11 A 55 + F 33 A^ 1 - 2F 31 A 45 ] 


F Z 
r 41 4 3 


45 


^ A 44 A 55 ” ^ A 45 * ^ ^ F 22 A 55 + F 44 A 44 " 2F 42 A 45 ^ 
E 0 = ^ A 44 A 55 ~ ^ A 45 ^ ^ ^- A 44 A 55 " ^ A 45 ^ ^ 
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where 


W 1 = F 33 ^ F 22 F 44 F 42 ^ F 32 F 44 + 2F 43 F 42 ' F 43 F 22 

W 2 = F 33 ^ F 22 A 55 + F 44 A 44 ” 2F 42 A 45 ^ ‘ A 55 ^ F 22 F 44 F 42 ^ 

+ F 32 A 55 “ 2F 43 F 32 A 45 + F 43 A 44 : 

2 

X 1 = F 31 ^ F 22 F 44 F 42 ^ F 32 ^ F 21 F 44 " F 41 F 42^ + F 43^ F 21 F 42 " F 41 F 22^ 


X 9 = -F 31 (F 22 A 55 ° + F aa A aa ° - 2F /19 A /ik 0 )- A,/ (F„F aa ~ F,/) 


44 44 


42 45 


1 / f F - F 2 ' 
-45 r 22 r 44 r 42 


+ F (F A °-F A °}-F (F A °-F A 

32 ' r 21 H 55 r 41 H 45 ' r 43 ir 21 H 45 r 41 H 44 ' 


Y 1 = F 31 (F 32 F 44 " F 43 F 42 ) “ F 33 (F 21 F 44 " F 41 F 42 ) + F 43 (F 21 F 43 “ F 41 F 32 ) 

Y 2 = '" A 45 1 (F 32 F 44 " F 43 F 42 ) " F 31 (F 32 A 55 " F 43 A 45 ) 

6 

1 / r- r- p“i~ v . ■ i— / r * n 0 


+ A 55 (F 21 F 44 “ F 41 F 42 ) + F 33 (F 21 A 55 


F 41 A 45 * 


Y 3 A 45* (F 32 A 55 “ F 43 A 45 0) " ^ (F 21 A 55 " F 41 A 45 ) 

Z 1 = F 31^ F 32 F 42 " F 43 F 22 ) “ F 33 (F 21 F 42 _ F 41 F 22 ) + F 32 (F 21 F 43 _ F 41 F 32^ 

Z 2 = F 31^ F 43 A 44 " F 32 A 45 ^ “ A 45 ^ F 32 F 42 “ F 43 F 22^ + A 55 ^ F 21 F 42 " F 41 F 22 

F 33^ F 41 A 44 " F 21 A 45 ^ 

Z 3 = " A 45 ^ F 43 A 44 " F 32 A 45 ^ + A 55 ^ F 41 A 44 ‘ F 21 A 45 ^ 


( 1 1 - 5 ) 
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when the coefficients A^, Bj g , Dj g (j=l,2) and A^ g are neglected 
we have 

v 13 = v 14 = U 11 = u 12 = 0 
and 

F 31 = F 45 = F 32 = F 42 = 0 ( 1 1-6) 

for this. case Equation (II -2 ) is 



and Equation (11) takes the uncoupled form 


( 1 1 “ 7 ) 
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By neglecting bending effects in sublaminates 1 and 0, the rotations 

p and p can be expressed in terms of p.. and p. and the characteristic 
oy ox ly xx 

equation reduces to 


(F llm s2 ‘ A 44 1) (F 22m s2 “ = 0 


where 


^11 h o 1 

F llm = T tA 22 (v ll + v 12 hT “T o * + B 22^ 


1 A 


44 


c h 2 rfl 1 , ^ h o A 55 . a 1, 

F 22m “ 2 [A 66 (u 13 u 14 h. . o B 66 ] 

1 a 55 


(II-9) 


Parameters F llm and F 22m represent a membrane behavior and therefore 


do not depend on the bending coefficients D-. 

«J 

Sublaminate 2: 


The elements of the matrix Cd in Equation (19) are defined as 


y 


[Cd] 


3x2 


A 


22 


'26 


22 




-1 




SYM 


A ^ 
H 12 

"z* 1 





B 66 1 

B 26 1 

°22 


B 12 1 



( 11 - 10 ) 


Displacements U 2 and V 2 and rotation p 2 ^ in sublaminate 2 are related 
to the applied uniform strain e and the rotation p 2x through the matrix Cd 
by 


43 



v 2,y 


e 

u 2,y 

= [Cd] 

P 2x,y 

_^2y,y 




( 11 - 11 ) 


Sublaminate 3: 


Parameters ^ 26 ancl ^66 ’ n Equation (24) are defined by 

J 22 = °22. + B 22 Wd 12 + B 26 wd 22 
J 26 = °26 + B 22 wd 13 + B 26 wd 23 
J 66 = D 66° + B 26° Wd 13 + B 66° wd 23 


where 


[wd] 


2x3 


A 0 A 0 

h 22 h 26 

A 0 A 0 

"26 H 66 


A 12° B 22 B 26 
A 16° B 26° B 66° 


(H-12) 

Displacements U 3 and V 3 in sublaminate 3 are related to the applied 


uniform strain e and the rotations 3 3x and p 3y through the matrix wd by 



[wd] 


P 3y,y 

P 3x,y 


(11-13) 
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Interlaminar Stresses: 


Parameters a., \p. and y. in Equation (29) are defined as 
J J J 


-1 


a . 
J 


(F 31 S j ” A 45 ) F 21 S j F 41 S j 


(F ll S j " A 44 ^ 

♦j 

= - 

(F 33 S j 2 ' A 55 1) F 32 S j' F 43 S j 


(F 31 S j 2 " A 45 ) 

_ r i_ 


* F 32 S j ^ F 22 S j “ A 44 * (F 42 S j " A 45 ]_ 


_ F 2* 2 J 


V j = V ll^ + a j V 13 + V 12 + Vl4 

U j = U 11 + a j U 13 + ^j U 12 + Y j U 14 

(11-14) 


The parameters associated with the resultant force and moment distributions in 
sublaminates 1 and 0, Equation (30) are given by 
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( 11 - 15 ) 


46 



Appendix III 


Classical Lamination Theory (CLT) Solution 
In the interior of the laminate CLT prevails that is, a membrane 
solution controls the behavior. The resultant forces N y and N are easily 
determined. In the following a ply-by-ply model is constructed. 

The constitutive relationship for each ply for a membrance behavior is 
given by. 



k 
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N x 


fl n 

A 12 

A 16 


£ 

N y 
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fl 12 

A 22 

A 26 


E y 

1 

X 

<< 

l 


Jl6 

A 26 

A 66 _ 


_ v. 


(III-l) 


Superscript k denotes the kth ply. According to the present formulation 




y x =U(y), y 


(III-2) 


Since e and y are functions of y only, therefore, from continuity of 

y x 

displacements at the interfaces between the kth and (k+l)th plies we have 


k+1 

k, h x k+1 , -h K \ 

u (y, 2 ) = u (y. — 2 — ) 


k, h k > k+1, -h k+1 , 

v (y, 2 ) = v (y. ~2 — ) 


(III-3) 


Substitute for the displacements from Equation (1) into Equation ( 1 1 1 -3 ) 
and differentiate w.r.t. y to get 
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(III-4) 


k _ k+1 

e y y 


• k _ k+i _ 

\ - Y x Y x 


From symmetry, consider one quarter of the laminate, 
of forces get 


From equilibrium 


N N . 

w k_V N k = 0 
M y " / xy 

k=l 


i » 

1 

k=l 


(III-5) 


Substitute from 


Equation (III-l) into (II1-5) to get 


S + $ 26 Y x + S 12 £ 


= 0 


22 y 

s 26 G y + S 66 Y x + S 16 " = ° 


where 


SLm 


N 

T A k 
L M> 

k=l 


H r 1,2 m = 2 and 6 


(1II-6) 


Solving for c y and y* from Equation (I1I-6) to get 

Oy * (s 26 s 16 " s 66 5 12 )e/S 
V* = < S 26 S 12 ‘ S 22 S 16 )e/S 

, a2 

s = s 22 S 66 " S 26 


(III ~7 ) 
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Substituting from Equation (III-7) into Equation ( 1 1 1 -1 ) to get 


N y k = (A 12 k e + A 22% + A 26 I<y x ) 

k k k k 

N = (A 1C e + A oc £ + A c , y ) 
xy . 16 26 y 66 'x' 


( I I 1-8) 


By comparing the expressions for C ancj C in Equation (33) with Equation 
(III -7 ) , we find that for sublaminates 1 and 0 


" C v = e y 


and 


C = y 
u 'x 


(III -9 ) 


Therefore the first term in Equation (30) for and in sublaminate 

1 and 0 is the CLT prediction. 

From equilibrium considerations the sign of the resultant axial force 

k 

Ny in each ply determines the sign of the interlaminar peel stress. This 
is illustrated for a [0/45/90/-451 with a delamination at the 0/45 inter- 
face. With reference to Figure 26, this layup corresponds to a 
[45/90/-45/0] laminate rotated at an angle 0 = -45 . 

T S 

* For the properties shown in Table I and using Equations ( 1 1 1 -7 ) and 
(III -8) the resultant forces N y and N xy in each ply can be determined. 
These are shown in Table XI. 

Consider the free body diagram of the top ply shown in Figure 29. From 
equilibrium of forces in the vertical direction the peel stress 
distribution should reverse its sign such that its resultant' is zero. 
Furthermore, the peel stress at the crack tip should be compressive in 
order to balance the moment of the compressive resultant force N^. 
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Table I. Material Properties and Geometry of the ED Specimen 


E 11 

= 134 

Gpa 

E 22 

= 10.2 

GPa 

G 12 

= 5.52 

GPa 

V 12 

= 0.3 


G 31 

= G 23 = 

3.1 


Ply thickness (h) = 0.14 x 10 3 m 

N 

Semi -Width (b) = 140h 
Crack length (a) = 6h 


Table II. Comparison of Characteristic Equation Coefficients 


Laminate 

-9 

Eg x 10 

i 

o 

X 

<X> 

UJ 

-9 

E 4 x 10 

-9 

E 2 x 10 

-9 

E x 10 
0 

[15^60/-75/-30] s 

0.143 

-6.940 x 10 3 

9.089 x 10 7 

-3.132 x 10 11 

3.193 x 10 14 

[-35/55^10/-80] s 

2.282 

-6.049 x 10 4 

4.614 x 10 8 

-1.155 x 10 12 

14 

5.676 x 10 1 

[40/-50/85(-5] s 

0.095 

-6.389 x 10 3 

1.071 x 10 8 

-4.096 x 10 11 

14 

3.193 x 10 1 


Table III. Comparison of the Characteristic Roots for a [35/80/-55/-101 Layup 

t s 


Model 



s 

x 2 

% 

Coupled 

5547 

4286 

2714 

2438 

Uncoupled 

5538 

3737 

2427 

2108 

Membrane 



2802 

2162 
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Table IV. COMPARISON OF ENERGY RELEASE COMPONENTS IN 
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Table V. [0/90/ + 45] g Laminate with Circular Hole- Edge Delamination at 
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Table Ml. [0/90/ + 45J s Laminate with Circular Hole- Edge Delamination at 
90/45 Interface (a = 6h, b = 140h. e = lO -6 ) 
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Table VII. [0/90/ + 45] Laminate with Circular Hole- Edge Delamination at 
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Table VIII. [45/90/-45/0] Laminate with Circular Hole- Edge Delamination at 
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Table IX. [45/90/-45/0] Laminate with Circular Hole- Edge Delamination at 
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Table X. [45/90/-45/0] Laminate with Circular Hole- Edge Delamination at 
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Table XI. Resultant Forces for a [0/45/90/-451 Laminate 

T S 


Ply Angle 

N y (N/m) 

N xy (N/m) 

0 

-0.008 

0 

45 

2.671 

3.031 

90 

-5.335 

0 

-45 

2.671 

-3.031 
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FIGURE 1. THE ED SPECIMEN 
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Fiber Orientation 
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FIGURE. 4 COMPARISON OF RELATIVE DISPLACEMENTS AT THE 

CRACK SURFACE FOR A [0/35/-35/90]s LAMINATE 
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FIGURE 5. COMPARISON OF RELATIVE DISPLACEMENTS AT THE 

CRACK SURFACE FOR A [-35/55/ I 0/-80] s LAMINATE 







FIGURE 7. NORMALIZED STRAIN-ENERGY-RELEASE 

RATE AS A FUNCTION OF DELAMINATION 
SIZE DIVIDED BY PLY THICKNESS 
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FIGURE 8. RESULTANT FORCE DISTRIBUTION 

IN A (30/-60/75/- 15) LAMINATE 
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FIGURE 9. RESULTANT FORCE DISTRIBUTION 

IN A (-35/55/10/80) LAMINATE 
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FIGURE 10. [0/90/+45] s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT +457-45 INTERFACE 
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FIGURE II. [0/90/±45] s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT +4S/-45 INTERFACE 
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FIGURE 12. [0/90/±45] s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT +457-45 INTERFACE 
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FIGURE 13. [0 
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FIGURE 14. [0/90/±45]s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT 90/+45 INTERFACE 
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FIGURE 15. [0/90/±45]s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT 90/+45 INTERFACE 



76 


£ 


FIGURE 16. [0/90/ ±45] s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT 90/+45 INTERFACE 
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FIGURE 17. [0/90/±45] s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT 0/90 INTERFACE 
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FIGURE 18. [0/90/±45]s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT 0/90 INTERFACE 
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FIGURE 19. [0/90/±45] s LAMINATE WITH CIRCULAR HOLE 

EDGE DELAMINATION AT 0/90 INTERFACE 
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FIGURE 20. COMPARISON OF ENERGY RELEASE RATE FOR 
[45/90/-45/0] s DELAMINATED AT -45/0 INTERFACE 
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FIGURE 21. COMPARISON OF ENERGY RELEASE RATE FOR 
[45/90/-45/0] s DELAMINATED AT -45/0 INTERFACE 
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FIGURE 22. COMPARISON OF ENERGY RELEASE RATE FOR 
[45/90/-45/0] s DELAMINATED AT -45/0 INTERFACE 
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FIGURE 26. COMPARISON OF ENERGY RELEASE RATE FOR 
[45/90/-45/0] s DELAMINATED -AT -45/0 INTERFACE 
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FIGURE 27. COMPARISON OF ENERGY RELEASE RATE FOR 
[45/90/-45/0] s DELAMINATED AT 45/90 INTERFACE 
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FIGURE 28. COMPARISON OF ENERGY RELEASE RATE FOR 
[45/90/-45/0] s DELAMINATED AT 45/90 INTERFACE 
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